.libPaths()
# delete data in global environment:
rm(list = ls())
# load standard packages:
library("stringr")
library("reshape")
library("ez")
library("ggplot2")
library("plyr")
library("nortest")
library("car")
library("ltm")
library("lsr")
library("afex")
library("nlme")
library("ggthemes")
library("robustHD")
library("cowplot")
library("tidyr")
library("lattice")                                
library("plyr")
library("DescTools")
library("data.table")

# load some necessary functions, for within-subject, repeated measure designs
load("/Users/XXX/effect.size.Rdata")
load("/Users/XXX/mean.sd.Rdata")
load("/Users/XXX/min.mean.sd.max.Rdata")
load("/Users/XXX/summarySE.Rdata")
load("/Users/XXX/summarySEwithin.Rdata")
load("/Users/XXX/normDataWithin.RData")
load("/Users/XXX/rcontrast.RData")

#read csv file:
getwd()
setwd("XXX")
data_for_r = read.csv("data_space_exp_future_banana.csv", fileEncoding="latin1")

#delete NAs
completeFun <- function(data, desiredCols) {
  completeVec <- complete.cases(data[, desiredCols])
  return(data[completeVec, ])
}
data_for_r<-completeFun(data_for_r, "rt")
# -> no NAs in the data

# only participants that have more than 50% correct answers per condition & answer:
data_for_r<-subset(data_for_r, 
                   participant==1 |
                     participant==2 |
                     participant==3 |
                     participant==8 |
                     participant==9 |
                     participant==11 |
                     participant==12 |
                     participant==15 |
                     participant==16 |
                     participant==17 |
                     participant==18 |
                     participant==19 |
                     participant==21 |
                     participant==23 |
                     participant==24 |
                     participant==25 |
                     participant==31 |
                     participant==33 |
                     participant==35 |
                     participant==37 |
                     participant==38 |
                     participant==42 |
                     participant==43 |
                     participant==45 |
                     participant==48 |
                     participant==49 |
                     participant==51 |
                     participant==53 |
                     participant==56 )

# response time: subtract 666ms (= stimulus length)
data_for_r$rt_offset<-data_for_r$rt-666
# only Responses that were given before stimulus onset (total amount of responses = 3480):
data_for_r<-subset(data_for_r, rt_offset>0) # 35 deleted
data_for_r<-subset(data_for_r, rt_offset<6000) # 54  deleted
# 89 (out of 3480) outliers = 2.5 %

#select only 'test blocks' (without learning blocks) per condition
data_for_r$test<-grepl(pattern = "(test)", data_for_r$blockcode) 
data_for_r<-subset(data_for_r, test=="TRUE")  

# rename conditions:
data_for_r<-mutate(data_for_r,blockcode=as.character(blockcode))
data_for_r<-mutate(data_for_r,blockcode=sapply(strsplit(data_for_r$blockcode, split='_', fixed=TRUE),function(x) (x[2])))
data_for_r$blockcode=gsub(data_for_r$blockcode, pattern = "vorwaerts", replacement = "forward")
data_for_r$blockcode=gsub(data_for_r$blockcode, pattern = "rueckwaerts", replacement = "backward")
data_for_r$blockcode=gsub(data_for_r$blockcode, pattern = "stehen", replacement = "standing")
data_for_r$condition <-data_for_r$blockcode

# only trials were an answer (correct or incorrect) was given
data_correct_and_incorrect<-subset(data_for_r, correct=="correct" | correct=="incorrect")  
# only correct answers:
data_correct<-subset(data_for_r, correct=="correct") 

#select relevant columns:
data_correct<-data_correct[c("participant", "Response", "condition", "correct", "rt_offset", "order", "trialnum", "blocknum", "MW_VVIQ2", "age")]

# count number of correct and incorrect items:
data_correct_count<-count(data_correct, c("correct", "participant", "condition", "Response", "order"))
data_correct_and_incorrect_count<-count(data_correct_and_incorrect, c("correct", "participant", "condition", "Response", "order"))

#### descriptive (summary for within-participant designs)
# per response:
dataR <- summarySEwithin(data_correct, measurevar="rt_offset", withinvars=c("Response"), idvar="participant")
dataR # -> rt of "ahead" items was around 100ms faster than "behind" items
# per condition:
dataCond <- summarySEwithin(data_correct, measurevar="rt_offset", withinvars=c("condition"), idvar="participant")
dataCond # -> no large differences between between conditions
# per response & condition
datac <- summarySEwithin(data_correct, measurevar="rt_offset", withinvars=c("condition", "Response"), idvar="participant")
datac
# -> "ahead" items during backward walking were around 200ms slower compared to "behind" items. 
# -> "ahead items during standing were around 50ms slower compared to "behind" items
# SDs from 425ms (ahead - forward) to 565ms (behind - backward)

## PLOT ##

ggplot(datac, aes(x=condition, y=rt_offset, fill=Response)) +
  geom_bar(position=position_dodge(.9), colour="black", stat="identity") +
  geom_errorbar(position=position_dodge(.9), width=.25, aes(ymin=rt_offset-ci, ymax=rt_offset+ci)) +
  theme_classic() +
  theme(legend.position="top") +
  xlab("\nCondition") +
  ylab("Response Time (ms)\n") +
  labs("d") +
  coord_cartesian(ylim=c(2500,1000)) +
  scale_fill_manual(values=c("#CCCCCC","#FFFFFF")) +
  theme_bw() +
  theme(text = element_text(size=30))+
  geom_hline(yintercept=38) 


## STATISTIK ##

# Check if response time correlates with trialnumber, blocknumber, age
cor(data_correct$rt_offset, data_correct$age)
plot(data_correct$rt_offset, data_correct$age)
cor(data_correct$rt_offset, data_correct$trialnum)
plot(data_correct$rt_offset, data_correct$trialnum)
cor(data_correct$rt_offset, data_correct$blocknum)
plot(data_correct$rt_offset, data_correct$blocknum)
cor(data_correct$rt_offset, data_correct$MW_VVIQ2)
# no correlations.

# Set contrast: Each moving condition compared to the control condition
# Therefore, the levels of 'condition' have to be re-ordered, with the control condition standing first. 
data_correct$condition <- factor(data_correct$condition, levels = c("standing", "backward", "forward"))
contrasts(data_correct$condition)

# Step 1: Baseline model 
m_baseline<-lme(rt_offset ~ 1, random = list(~1|participant, ~1+condition|participant, ~1+Response|participant), method = "ML", data = data_correct)

# Step 2: condition
m_condition<-update(m_baseline, ~ + condition)
anova(m_baseline, m_condition) # -> n.s.

# Step 3: Response
m_Response<-update(m_baseline, ~ + Response)
anova(m_baseline, m_Response) # significant

# Step 4: condition + Response
m_cond_and_res<-update(m_baseline, ~ + Response + condition)

# Step 5: Interaction condition * Response
m_cond_and_res_interaction<-update(m_baseline, ~ + Response * condition)
anova(m_cond_and_res,m_cond_and_res_interaction) # -> significant
summary(m_cond_and_res_interaction)

#Check assumptions:
qqnorm(m_cond_and_res_interaction)
plot(m_cond_and_res_interaction)
# -> ok 

# lme model adding order
lme_model_Order<-update(m_cond_and_res_interaction, .~. + order)
summary(lme_model_Order)
anova(m_cond_and_res_interaction, lme_model_Order) # -> n.s.

# post hoc tests 
data_backward<-subset(data_correct, data_correct$condition =="backward")
backward_post<-aggregate(rt_offset ~ Response + participant, data=data_backward, FUN=mean)
backward_post$participant<-as.factor(backward_post$participant)
pairedSamplesTTest(rt_offset ~ Response, data=backward_post, id = "participant")
# -> p = .01, Cohen's d:.49

data_forward<-subset(data_correct, data_correct$condition =="forward")
forward_post<-aggregate(rt_offset ~ Response + participant, data=data_forward, FUN=mean)
forward_post$participant<-as.factor(forward_post$participant)
pairedSamplesTTest(rt_offset ~ Response, data=forward_post, id = "participant")
# -> no differences

data_standing<-subset(data_correct, data_correct$condition =="standing")
standing_post<-aggregate(rt_offset ~ Response + participant, data=data_standing, FUN=mean)
standing_post$participant<-as.factor(standing_post$participant)
pairedSamplesTTest(rt_offset ~ Response, data=standing_post, id = "participant")
# -> no differences

